#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
#@@   CHOOSE THE INPUT VARIABLES       @@#
#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
compile="no"
DO_VALENCE_ANALISIS="yes"
	Do_Zq_rprom_and_qav="yes"
        Do_D_rprom_and_qav="yes"
        Do_ChiralVal="yes"
#        Do_prom_pm_chiral="yes"
#                # Optional:
#                Do_RGI_valence="yes"
#                Do_p2_val_fit="yes"
DO_CHIRAL_SEA="no"
DO_RGI="no"
DO_P2FIT="no"


source ~/workdir/root/bin/thisroot.sh

##########################
#Beta (beta*100)
#beta100=190
#beta100=195
#beta100=210
beta100=${1}
beta=$(echo "scale=2; $beta100/100" | bc -l)


#############################

################################
## Ensemble List
if [ $beta100 -eq 195 ]
then
#eList=(1 2 3 4 7 8 )
eList=(1 2 3 7 8)
fi
if [ $beta100 -eq 190 ]
then
#eList=("A1")
eList=("A3")
fi
if [ $beta100 -eq 210 ]
then
eList=(2 3 4 5)
fi
#####################################
#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
#   For CHIRAL SEA analysis steps choose e=0
#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
eList=("${2}")

#Choose theta: p, m
#pmList=("p" "m")
#pmList=("m")
#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
#   For steps after p/m average choose pm=a
#!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
pmList=("${3}")
###################################

#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
#@@   FIXED INPUT VARIABLES    
#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#


###################################
# Directories
dirreno="/home/nucave/workdir/lattice/RENORMALIZATION/reno"

dir="$dirreno/NF4/Z4F/"
dirinput="$dirreno/NF4/Z4F/INPUT"

dirdata="/home/nucave/workdir/lattice/RENORMALIZATION/"
######################################

# Nclust
Nclust=16
#Gauge Fixing and Gauge Action
gauge_fixing="landau"
gauge_action="iwasaki"

# Fermionic Action
Nf=4
action="TwistedMass"

#Nr=2 -> r0 & r1
Nr=2

#@@   depending on beta          
#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#

#####################################
#Gamma and Momentum order convention
#WARNIG: No correr a la vez scripts para diferentes definiciones de gammas o momentos

if [ $beta100 -eq 190 ] || [ $beta100 -eq 210 ]
then
gamma="TWISTED"
momentum_order="TYZX"
convention="TW"
fi

if [ $beta100 -eq 195 ]
then
gamma="APE"
momentum_order="TYZX"
convention="APE"
fi

# Original Francesco format
#gamma="TWISTED"
#momentum_order="TXYZ"
#convention="HUYGENS-CINECA"

# L^3xT=V
if [ $beta100 -eq 190 ]
then
 L=24
 T=48
fi
if [ $beta100 -eq 195 ]
then
 L=24
 T=48
fi
if [ $beta100 -eq 210 ]
then
 L=32
 T=64
fi

##################################
# Oa2 corrections
Oa2List=("uncor" "Oa2g2" "Oa2g2b")

NOa2=${#Oa2List[*]}

if [ $beta100 -eq 190 ]; then  Uplaq=0.50; fi
if [ $beta100 -eq 195 ]; then  Uplaq=0.58; fi
if [ $beta100 -eq 210 ]; then  Uplaq=0.6133375; fi

##########################################

######################################3
# ainv for RGI
if [ $beta100 -eq 190 ]; then  ainv=2.19 ; fi
if [ $beta100 -eq 195 ]; then  ainv=2.50 ; fi
if [ $beta100 -eq 210 ]; then  ainv=3.23 ; fi

############## This will be the common case

#if [ $Do_RGI_valence = "no" ]
#then
#rgiList=(\"\")
#fi

########## Just for intermediate checks:
#if [ $Do_RGI_valence = "yes" ]
#then
#rgiList=(\"\")

# If I have computed RGI also on the valence I can read directly these rgi files and don't do the RGI on the sea.
#rgiList=("_RIp_ainv")

#fi



##################################
## p2 FIT

# p2 cuts
cut=0.28
#cut=100

p2min_lin=1.2
p2max_lin=2.2

if [ $beta100 -eq 195 ]
then
p2min_av=1.8
p2max_av=2.0
fi

if [ $beta100 -eq 210 ]
then
p2min_av=1.1
p2max_av=1.23
fi
if [ $beta100 -eq 190 ]
then
p2min_av=2.4
p2max_av=2.7
fi
######################################

#####################################
# Momenta
## $dirinput/mom file

if [ $beta100 -eq 190 ]
then
Nmom=492
Nmom_eq=135
fi

if [ $beta100 -eq 195 ]
then
Nmom=140
Nmom_eq=55

fi
if [ $beta100 -eq 210 ]
then
Nmom=492
Nmom_eq=135
fi


echo "1 0 0 0" > $dirinput/mom_${beta100}.dat # boundary conditions
echo "2" >> $dirinput/mom_${beta100}.dat # Nindexmom
if [ $Nmom -eq 140 ]
then
echo "0 3 0 2 0 2 0 2" >> $dirinput/mom_${beta100}.dat
echo "4 7 2 3 2 3 2 3" >> $dirinput/mom_${beta100}.dat
fi
if [ $Nmom -eq 492 ]
then
echo "0 3 0 2 0 2 0 2" >> $dirinput/mom_${beta100}.dat
echo "4 9 2 5 2 5 2 5" >> $dirinput/mom_${beta100}.dat
fi

#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#

if [ $compile = "yes" ]
then
echo "Compiling..."

###### If you're going to launch a script don't compile. 
###### You may have different scripts compiling the same program
###### Compile before or do it just the first time.

sh compile.sh D_rPROM_QAVERAGE
sh compile.sh Subst_Oa2g2_Z4F
sh compile.sh D_pmAVERAGE
sh compile.sh Calc_Z_oneP
sh compile.sh Z4FRGI

fi

if [ $DO_VALENCE_ANALISIS = "yes" ]
then


#@@@@@@@@@@@ BEGIN ENSEMBLE LIST

for e in  ${eList[*]}
do
                echo $e
                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
                #@@  MASSES TO BE ANALYSED  @@#
                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#

                if [ $beta100 -eq 190 ]
                then
                        if [ $e = "A1" ]; then Nmass=6; fi
                        if [ $e = "A3" ]; then Nmass=6; fi
                        if [ $e = "A4" ]; then Nmass=6; fi
                fi


                if [ $beta100 -eq 195 ]
                then
                        if [ $e -eq 1 ] || [ $e -eq 2 ] || [ $e -eq 7 ] || [ $e -eq 8 ] ; then Nmass=5 ; fi
                        if [ $e -eq 3 ] || [ $e -eq 4 ] ; then Nmass=8 ; fi
                fi

                if [ $beta100 -eq 210 ]
                then
                        if [ $e -eq 2 ]; then Nmass=7; fi
                        if [ $e -eq 3 ]; then Nmass=7; fi
                        if [ $e -eq 5 ]; then Nmass=6; fi
                        if [ $e -eq 4 ]; then Nmass=6; fi
                fi
                Nmassm1=$( echo "$Nmass -1"| bc -l)


        for pm in ${pmList[*]}
        do
                epm=${e}${pm}
                echo $epm 
                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
                #@@  Configurations to be read:
                #@@  conf=conf_ini+i*jump with 
                #@@      i=0,1,...Nconf-1
                #@@      Nconf=Nex*Nclust
                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#


                if [ $beta100 -eq 190 ]
                then

                        if [ $epm = "A1m" ] ; then conf_ini=129;  jump=6; Nex=15; fi
                        if [ $epm = "A1p" ] ; then conf_ini=499 ; jump=6; Nex=15; fi

                        if [ $epm = "A3m" ] ; then conf_ini=198 ; jump=6 ; Nex=15; fi
                        if [ $epm = "A3p" ] ; then conf_ini=499 ; jump=6 ; Nex=15; fi

                        if [ $epm = "A4m" ] ; then conf_ini=499 ; jump=6; Nex=15; fi
                        if [ $epm = "A4p" ] ; then conf_ini=499 ; jump=4; Nex=20; fi


                fi


                if [ $beta100 -eq 195 ]
                then
                        jump=4

                        if [ $e -eq 1 ] ; then  Nex=19 ;fi
                        if [ $e -eq 2 ] ; then  Nex=22 ;fi
                        if [ $e -eq 3 ] ; then  Nex=22 ;fi
                        if [ $e -eq 4 ] ; then  Nex=14 ;fi
                        if [ $e -eq 7 ] ; then  Nex=23 ;fi
                        if [ $e -eq 8 ] ; then  Nex=19 ;fi

                        if [ $epm = "1m" ] ; then conf_ini=118 ;fi
                        if [ $epm = "1p" ] ; then conf_ini=102 ;fi

                        if [ $epm = "2m" ] ; then conf_ini=105 ;fi
                        if [ $epm = "2p" ] ; then conf_ini=384 ;fi

                        if [ $epm = "3m" ] ; then conf_ini=100 ;fi
                        if [ $epm = "3p" ] ; then conf_ini=241 ;fi

# WARNING WARNING e=4 There are corrupted confs!
                        if [ $epm = "4m" ] ; then conf_ini=459 ;fi
                        if [ $epm = "4p" ] ; then conf_ini=172 ;fi

                        if [ $epm = "7m" ] ; then conf_ini=1454 ;fi
                        if [ $epm = "7p" ] ; then conf_ini=102 ;fi

                        if [ $epm = "8m" ] ; then conf_ini=127 ;fi
                        if [ $epm = "8p" ] ; then conf_ini=326 ;fi
                fi
                if [ $beta100 -eq 210 ]
                then

                        if [ $epm = "2m" ] ; then jump=4; Nex=12; conf_ini=260; fi

                        if [ $epm = "3p" ] ; then jump=1; Nex=15; conf_ini=100; fi
                        if [ $epm = "3m" ] ; then jump=1; Nex=22; conf_ini=100; fi
                        if [ $epm = "4p" ] ; then jump=1; Nex=13; conf_ini=100; fi
                        if [ $epm = "4m" ] ; then jump=1; Nex=23; conf_ini=100; fi
                        if [ $epm = "5p" ] ; then jump=1; Nex=14; conf_ini=100; fi
                        if [ $epm = "5m" ] ; then jump=1; Nex=18; conf_ini=100; fi

                fi


                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
                #@@  kappa_sea                       
                #@   (only for name construction)
                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#

                if [ $beta100 -eq 195 ]
                then
                        if [ $epm = "1m" ] ; then kappa_sea="0.160754" ;fi
                        if [ $epm = "1p" ] ; then kappa_sea="0.162145" ;fi

                        if [ $epm = "2m" ] ; then kappa_sea="0.161229" ;fi
                        if [ $epm = "2p" ] ; then kappa_sea="0.160826" ;fi

                        if [ $epm = "3m" ] ; then kappa_sea="0.161229" ;fi
                        if [ $epm = "3p" ] ; then kappa_sea="0.160826" ;fi

                        if [ $epm = "4m" ] ; then kappa_sea="0.161095" ;fi
                        if [ $epm = "4p" ] ; then kappa_sea="0.16087" ;fi

                        if [ $epm = "7m" ] ; then kappa_sea="0.160681" ;fi
                        if [ $epm = "7p" ] ; then kappa_sea="0.161925" ;fi

                        if [ $epm = "8m" ] ; then kappa_sea="0.160681" ;fi
                        if [ $epm = "8p" ] ; then kappa_sea="0.161925" ;fi

                fi
                if [ $beta100 -eq 210 ]
                then
                        if [ ${epm} = "2m" ] ; then kappa_sea="0.156157" ; fi

                        if [ ${epm} = "3m" ] ; then kappa_sea="0.156209" ; fi
                        if [ ${epm} = "3p" ] ; then kappa_sea="0.156017" ; fi
                        if [ ${epm} = "4m" ] ; then kappa_sea="0.156250" ; fi
                        if [ ${epm} = "4p" ] ; then kappa_sea="0.155983" ; fi
                        if [ ${epm} = "5m" ] ; then kappa_sea="0.156291" ; fi
                        if [ ${epm} = "5p" ] ; then kappa_sea="0.155949" ; fi
                fi
                if [ $beta100 -eq 190 ]
                then
                        if [ ${epm} = "A1p" ] ; then kappa_sea="0.163206" ; fi

                        if [ ${epm} = "A3p" ] ; then kappa_sea="0.164537" ; fi
                        if [ ${epm} = "A3m" ] ; then kappa_sea="0.162260" ; fi
                        if [ ${epm} = "A4p" ] ; then kappa_sea="0.164817" ; fi
                        if [ ${epm} = "A4m" ] ; then kappa_sea="0.162113" ; fi

                fi


                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
                #@@ INPUT OUTPUT DIRECTORIES
                #@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
		
		if [ $convention = "APE" ]
		then
        	      dir_data="$dirdata/DATA/OUTPUT-APE/${action}Nf${Nf}/b$beta100/Z4F/${e}${pm}/"
                      dir_data_doublebarra="\/home\/nucave\/workdir\/lattice\/RENORMALIZATION\/DATA\/OUTPUT-APE\/${action}Nf${Nf}\/b$beta100\/Z4F\/${e}${pm}\/"

		else
		      dir_data="$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/${e}${pm}/"
                      dir_data_doublebarra="\/home\/nucave\/workdir\/lattice\/RENORMALIZATION\/DATA\/OUTPUT\/${action}Nf${Nf}\/b$beta100\/Z4F\/${e}${pm}\/"

		fi
		mkdir -p -v $dir_data/QAV/Z_ALLP
		mkdir -p -v $dir_data/QAV/D_ALLP
		for i in `seq 0 $Nmom_eq`; do imom=`printf %03d $i`; mkdir -p -v $dir_data/QAV/p$imom ; done

		cd $dir/EXECUTABLES/
		#######################

			# 1.) Zq r-prom and q2-av
		 if [ $Do_Zq_rprom_and_qav = "yes" ]
		 then
			time ./Subst_Oa2g2_rPROM_QAVERAGE_Zq -L=$L -T=$T  -dir=$dir_data -kappa_sea=$kappa_sea -conf_ini=${conf_ini} -jump=${jump} -Nclust=${Nclust} -Nex=${Nex}  -Nr=${Nr} -beta100=${beta100}  -NOa2=${NOa2}  -Uplaq=${Uplaq} -beta=${beta}  -gauge_fixing=${gauge_fixing} -gauge_action=${gauge_action} -Nf=${Nf} -Nmass=${Nmass} -momentum_order=${momentum_order}
		fi

			# 2.) D r-prom and q2-av
		if [ $Do_D_rprom_and_qav = "yes" ]
		then
		       time  for im in `seq 0 $Nmassm1`; do ./D_rPROM_QAVERAGE -L=$L -T=$T  -dir=$dir_data -kappa_sea=$kappa_sea -conf_ini=${conf_ini} -jump=${jump} -Nclust=${Nclust} -Nex=${Nex}  -Nr=${Nr} -beta100=${beta100} -imass=${im} -momentum_order=${momentum_order}; done
		fi

			# 3.) Chiral Valence limit

	if [ $Do_ChiralVal = "yes" ]
	then

			# Masses values
                epm=${e}${pm}
                if [ $beta100 -eq 190 ]
                then
                        if [ $e = "A1" ] || [ $e = "A3" ] || [ $e = "A4" ]
                        then
                                muList=(0.0060 0.0080 0.0120 0.0170 0.0210 0.0260)
                        fi

                        if [ $epm = "A1p" ];then mPCAC=0.01629; fi
                        if [ $epm = "A1m" ];then mPCAC=0.01588; fi
                        if [ $epm = "A1a" ];then mPCAC=0.01085; fi

                        if [ $epm = "A2p" ];then mPCAC=0.01552; fi
                        if [ $epm = "A2m" ];then mPCAC=0.01536; fi
                        if [ $epm = "A2a" ];then mPCAC=0.01544; fi

                        if [ $epm = "A3p" ];then mPCAC=0.01525; fi
                        if [ $epm = "A3m" ];then mPCAC=0.01469; fi
                        if [ $epm = "A3a" ];then mPCAC=0.01497; fi

                        if [ $epm = "A4p" ];then mPCAC=0.0142; fi
                        if [ $epm = "A4m" ];then mPCAC=0.0147; fi
                        if [ $epm = "A4a" ];then mPCAC=0.01445; fi


                        if [ $epm = "A1p" ]; then MPSList=(0.268120752 0.273678855 0.286937707 0.305098627 0.319859691 0.338077073 0.278840807 0.291457274 0.309043702 0.32347193 0.341377957 0.302972885 0.319507306 0.333291525 0.350567241 0.334955306 0.348031935 0.364580305 0.360625206 0.376658069 0.392125341 );fi
                        if [ $epm = "A1m" ]; then MPSList=(0.268120752 0.273678855 0.286937707 0.305098627 0.319859691 0.338077073 0.278840807 0.291457274 0.309043702 0.32347193 0.341377957 0.302972885 0.319507306 0.333291525 0.350567241 0.334955306 0.348031935 0.364580305 0.360625206 0.376658069 0.392125341 );fi


                        if [ $epm = "A2p" ]; then MPSList=(0.267229315 0.273170904 0.287141557 0.305988268 0.321164316 0.339785076 0.278659948 0.291911167 0.310124191 0.324935897 0.34321821 0.303947782 0.321014045 0.335128372 0.35273204 0.336909039 0.350266283 0.367097069 0.36311123 0.379400286 0.395100708); fi
                        if [ $epm = "A2m" ]; then MPSList=(0.224203831 0.231219298 0.246401662 0.265887977 0.281275929 0.300013133 0.237639537 0.251881371 0.270508154 0.285370339 0.3035813 0.264568356 0.281720803 0.295663247 0.312944347 0.29741693 0.310423501 0.326738573 0.322812903 0.338477587 0.353444678) ; fi

                        if [ $epm = "A3p" ]; then MPSList=(0.265216934 0.271531587  0.286146332 0.30558633 0.321123553  0.340114756 0.277357728  0.291219721 0.31002186 0.325202175 0.343866774 0.303832602 0.321482491 0.335973696 0.353973709 0.337964303 0.351704869 0.368939947 0.364936478 0.381631347 0.397734775); fi
                        if [ $epm = "A3m" ]; then MPSList=( 0.217781488 0.224794155  0.240061987  0.259597882 0.274955528 0.293588959 0.231201294 0.245534713 0.264234336 0.279092049 0.297240315 0.258339699 0.275603882 0.289577184 0.306845052 0.291453023 0.30452239 0.320863707 0.316996687 0.332712234 0.347755052) ; fi

                        if [ $epm = "A4p" ]; then MPSList=(0.267607195  0.274105196 0.289091476  0.308918239 0.324696294  0.34391784  0.280074906  0.294250278 0.313389665 0.328783706 0.347656249 0.30709032 0.325002607 0.339664686 0.35783446 0.341680667 0.35555297 0.372920968 0.368892498 0.385697234 0.401887928);  fi
                        if [ $epm = "A4m" ]; then MPSList=(0.212383364 0.220162098 0.236317373 0.256514771 0.272275847 0.291343384 0.227169963 0.242205796 0.261450014 0.27664595 0.295153440 0.255502533 0.273170158 0.287398407 0.304931356 0.289317522 0.302575237 0.319106853 0.315195311 0.331055742 0.346198771 ); fi

                fi

                if [ $beta100 -eq 195 ]
                then
                        if [ $e -eq 1 ] || [ $e -eq 2 ] || [ $e -eq 7 ] || [ $e -eq 8 ] ; then muList=(0.0085 0.015 0.0203 0.0252 0.0298) ; fi
                        if [ $e -eq 3 ] || [ $e -eq 4 ] ; then muList=(0.0060 0.0085 0.0120 0.0150 0.0180 0.0203 0.0252 0.0298) ; fi

                        if [ $epm = "1p" ];then mPCAC=0.0195; fi # * mPCAC_val ! = mPCAC_sea
                        if [ $epm = "1m" ];then mPCAC=-0.0216; fi
                        if [ $epm = "1a" ];then mPCAC=0.0200; fi

                        if [ $epm = "2p" ];then mPCAC=0.0191; fi
                        if [ $epm = "2m" ];then mPCAC=-0.0213; fi
                        if [ $epm = "2a" ];then mPCAC=0.0200; fi

                        if [ $epm = "3p" ];then mPCAC=0.0162; fi
                        if [ $epm = "3m" ];then mPCAC=-0.0160; fi
                        if [ $epm = "3a" ];then mPCAC=0.0161; fi

                        if [ $epm = "4p" ];then mPCAC=0.0151; fi
                        if [ $epm = "4m" ];then mPCAC=-0.0146; fi
                        if [ $epm = "4a" ];then mPCAC=0.0149; fi

                        if [ $epm = "7p" ];then mPCAC=0.0181; fi # * mPCAC_val ! = mPCAC_sea
                        if [ $epm = "7m" ];then mPCAC=-0.0180; fi
                        if [ $epm = "7a" ];then mPCAC=0.0180; fi

                        if [ $epm = "8p" ];then mPCAC=0.0183; fi # * mPCAC_val ! = mPCAC_sea
                        if [ $epm = "8m" ];then mPCAC=-0.0194; fi
                        if [ $epm = "8a" ];then mPCAC=0.0187; fi




                        if [ $epm = "1p" ] ; then MPSList=(0.2955 0.3127 0.3290 0.3445 0.3591 0.3282 0.3434 0.3581 0.3721 0.3578 0.3719 0.3855 0.3856 0.3988 0.4116) ;fi
                        if [ $epm = "1m" ] ; then MPSList=(0.2553 0.2732 0.2895 0.3050 0.3195 0.2888 0.3037 0.3181 0.3318 0.3177 0.3313 0.3443 0.3444 0.3569 0.3690) ;fi

                        if [ $epm = "2p" ] ; then MPSList=(0.2764 0.2925 0.3079 0.3226 0.3365 0.3072 0.3216 0.3357 0.3491 0.3354 0.3489 0.3619 0.3621 0.3747 0.3870) ;fi
                        if [ $epm = "2m" ] ; then MPSList=(0.2587 0.2770 0.2937 0.3094 0.3240 0.2930 0.3082 0.3227 0.3366 0.3224 0.3362 0.3494 0.3494 0.3622 0.3745) ;fi

                        if [ $epm = "3p" ] ; then MPSList=(0.2550 0.2605 0.2701 0.2792 0.2886 0.2960 0.3117 0.3263 0.2658 0.2749 0.2836 0.2928 0.3000 0.3154 0.3298 0.2834 0.2918 0.3006 0.3076 0.3226 0.3366 0.2998 0.3084 0.3151 0.3298 0.3436 0.3167 0.3233 0.3376 0.3512 0.3297 0.3438 0.3572 0.3575 0.3705 0.3832) ;fi
                        if [ $epm = "3m" ] ; then MPSList=(0.1931 0.2062 0.2233 0.2370 0.2500 0.2597 0.2792 0.2966 0.2179 0.2335 0.2464 0.2587 0.2678 0.2865 0.3033 0.2477 0.2594 0.2709 0.2795 0.2973 0.3134 0.2705 0.2814 0.2896 0.3068 0.3223 0.2918 0.2997 0.3163 0.3314 0.3074 0.3235 0.3383 0.3389 0.3531 0.3668) ;fi

                        if [ $epm = "4p" ] ; then MPSList=(0.2434 0.2495 0.2598 0.2694 0.2792 0.2869 0.3031 0.3181 0.2552 0.26550 0.2742 0.2838 0.2913 0.3073 0.3220 0.2742 0.2830 0.2922 0.2994 0.3149 0.3294 0.2915 0.3004 0.3074 0.3226 0.3368 0.3091 0.3160 0.3307 0.3446 0.3227 0.3372 0.3509 0.3513 0.3646 0.3775) ;fi
                        if [ $epm = "4m" ] ; then MPSList=(0.2140 0.2229 0.2365 0.2483 0.2599 0.2687 0.2869 0.3034 0.2309 0.2433 0.2543 0.2653 0.2737 0.2912 0.3072 0.2545 0.2646 0.2749 0.2828 0.2994 0.3148 0.2741 0.2839 0.2914 0.3075 0.3224 0.2932 0.3005 0.3161 0.3305 0.3076 0.3228 0.3369 0.3372 0.3509 0.3640) ;fi

                        if [ $epm = "7p" ] ; then MPSList=(0.282638091 0.301027689 0.317993738 0.334000974 0.348939203 0.317524386 0.333327372 0.348493733 0.362795744 0.348338492 0.362892201 0.37670548  0.376950823 0.390356011 0.403416707) ;fi
                        if [ $epm = "7m" ] ; then MPSList=( 0.238559831  0.259201635 0.277090595 0.293612613 0.308888896 0.27693825  0.293093636 0.308342572 0.322625194 0.308164817  0.3225951  0.336226335  0.336393132  0.349512089  0.362209606) ;fi

                         if [ $epm = "8p" ] ; then MPSList=(0.282924605  0.301179812  0.318035273  0.333948003  0.348803474  0.317633898  0.333377744  0.348483981  0.362727689  0.348357462  0.362869163  0.376637287  0.376897984  0.390267435  0.403300018) ;fi
                         if [ $epm = "8m" ] ; then MPSList=(0.244623233  0.264032004  0.281343342  0.297516075  0.312557041  0.280883845  0.296597472  0.311564494  0.32565134   0.311263273  0.325426169  0.338871341  0.338964708  0.351904801  0.364429301) ;fi


                fi

                if [ $beta100 -eq 210 ]
                then

                        if [ $e -eq 2 ]
                        then
                                muList=(0.0013 0.0030 0.0080 0.0143 0.0195 0.0247 0.0298)

                        fi

                        if [ $e -eq 3 ]
                        then
                                muList=(0.0025 0.0046 0.0090 0.0152 0.0201 0.0249 0.0297)

                        fi
                        if [ $e -eq 5 ]
                        then
                                muList=(0.0048  0.0078  0.0119  0.0190  0.0242  0.0293)
                        fi
                        if [ $e -eq 4 ]
                        then
                                muList=(0.0039 0.0064 0.0184 0.0240 0.0295)
                        fi


                        # CHECK en la tabla del POS no esta ni 2m/p ni 4m/p

                       if [ ${epm} = "2a" ];then mPCAC=0.0041; fi
                       if [ ${epm} = "2m" ];then mPCAC=-0.00403; fi
                       if [ ${epm} = "2p" ];then mPCAC=0.0042; fi

                       if [ ${epm} = "3a" ];then mPCAC=0.00572; fi
                       if [ ${epm} = "3m" ];then mPCAC=-0.00585; fi
                       if [ ${epm} = "3p" ];then mPCAC=0.00559; fi

                       if [ ${epm} = "4a" ];then mPCAC=0.00685; fi
                       if [ ${epm} = "4m" ];then mPCAC=-0.00682; fi
                       if [ ${epm} = "4p" ];then mPCAC=0.00684; fi

                       if [ ${epm} = "5a" ];then mPCAC=0.0082; fi
                       if [ ${epm} = "5m" ];then mPCAC=-0.0082; fi
                       if [ ${epm} = "5p" ];then mPCAC=0.0082; fi


                fi

                ZAPREV=0.77
                Nmass=${#muList[*]}
                Nmassm1=$(echo "$Nmass - 1 " | bc -l )
                Nmasscomb=${#MPSList[*]}

                #################################
                ## INPUT OF THE PROGRAMS
                ################################
                # mu_$epm.dat
                rm -f  $dirinput/mu_${e}${pm}_${beta100}.dat
                echo $Nmass > $dirinput/mu_${e}${pm}_${beta100}.dat
                for im in `seq 1 $Nmass`
                do
                                mass=${muList[$(echo " $im -1 " |bc -l )]}
                                echo $mass >> $dirinput/mu_${e}${pm}_${beta100}.dat
                done
                echo " ---> $dirinput/mu_${e}${pm}_${beta100}.dat file written"

                rm -f  $dirinput/MPS_${e}${pm}_${beta100}.dat
                for im in `seq 1 $Nmasscomb`
                do
                                MPS=${MPSList[$(echo " $im -1 " |bc -l )]}
                                echo $MPS >> $dirinput/MPS_${e}${pm}_${beta100}.dat
                done
                echo "  --> $dirinput/MPS_${e}${pm}_${beta100}.dat file written"        

                # mPCAC_$epm.dat
                rm -f $dirinput/mPCAC_${e}${pm}_${beta100}.dat
                echo $mPCAC > $dirinput/mPCAC_${e}${pm}_${beta100}.dat
                echo $ZAPREV >> $dirinput/mPCAC_${e}${pm}_${beta100}.dat
                echo "  --> $dirinput/mPCAC_${e}${pm}_${beta100}.dat  file written"        

	        if [ $convention = "APE" ]
                then

                dir_output_doblebarra="\/home\/nucave\/workdir\/lattice\/RENORMALIZATION\/DATA\/OUTPUT-APE\/${action}Nf${Nf}\/b$beta100\/Z4F\/${e}${pm}\/"
                else
                dir_output_doblebarra="\/home\/nucave\/workdir\/lattice\/RENORMALIZATION\/DATA\/OUTPUT\/${action}Nf${Nf}\/b$beta100\/Z4F\/${e}${pm}\/"
                fi
	
		#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
		#@@@ CHIRAL LIMIT FOR Zq       @@@@@@@
		#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
	        for Oa2print in ${Oa2List[*]}
        	do

	        time for imom in `seq 0 $Nmom_eq`
                do
                        echo "=================================================================="
                        echo "=======>         $Oa2print p=$imom             <=================="

                                sed "s/__MOMENTUM_ORDER__/$momentum_order/g;s/__BETA100__/$beta100/g;s/__L__/$L/g;s/__T__/$T/g;s/__NCLUST__/$Nclust/g;s/__DIROUTPUT__/$dir_output_doblebarra/g;s/__IMOM__/$imom/g;s/__ENSEMBLE__/${e}/g;s/__PMI__/${pm}/g;s/__Oa2PRINT__/${Oa2print}/g;s/__GAMMA__/"q"/g;s/__IPROP__/${iprop}/g;;s/__USEMPS__/0/g;s/__SUBGP__/0/g" ../PROGRAMS/Z_CHIRALVAL_template > ../EXECUTABLES/Z_CHIRALVAL_${beta100}_${e}${pm}
                                root -q -b -l ../EXECUTABLES/Z_CHIRALVAL_${beta100}_${e}${pm}
                                rm ../EXECUTABLES/Z_CHIRALVAL_${beta100}_${e}${pm}
                done
		done
		
		#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@#
		#@@@@ CHIRAL LIMIT FOR D                        @@@@
		#@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@@
	        time for imom in `seq 0 $Nmom_eq`
		do
	
			  # GP2 = 1 -> Dij with i>=2 is fit with a double GP
			  # b=1.95 -> Nmass=5 Si hago el polo doble son 4 parametros en el fit por lo que es dificil ver el polo doble 
			  # Hago polo simple
			  # USEMPS=1 -> Substract GP with MPS^2 instead of mpolar

			  sed "s/__MOMENTUM_ORDER__/$momentum_order/g;s/__BETA100__/$beta100/g;s/__L__/$L/g;s/__T__/$T/g;s/__NCLUST__/$Nclust/g;s/__DIR__/$dir_data_doublebarra/g;s/__IMOM__/$imom/g;s/__ENSEMBLE__/${e}/g;s/__PM__/${pm}/g;s/__USEMPS__/1/g;s/__GP2__/1/g" ../PROGRAMS/Z4F_CHIRALVAL_template > ../EXECUTABLES/Z4F_CHIRALVAL_${beta100}_${e}${pm} 
			 root -q -b -l ../EXECUTABLES/Z4F_CHIRALVAL_${beta100}_${e}${pm}
			 rm ../EXECUTABLES/Z4F_CHIRALVAL_${beta100}_${e}${pm}  
		 done

		fi

	done # pmList




	#######################
	## 4.) pm-average
	########################
        if [ $convention = "APE" ]
        then
              dir_av="$dirdata/DATA/OUTPUT-APE/${action}Nf${Nf}/b$beta100/Z4F/${e}a/"
	      dir_plus="$dirdata/DATA/OUTPUT-APE/${action}Nf${Nf}/b$beta100/Z4F/${e}p/"
	      dir_minus="$dirdata/DATA/OUTPUT-APE/${action}Nf${Nf}/b$beta100/Z4F/${e}m/"
        else
              dir_av="$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/${e}a/"
	      dir_plus="$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/${e}p/"
	      dir_minus="$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/${e}m/"
        fi

        mkdir -p -v $dir_av/QAV/Z_ALLP
        mkdir -p -v $dir_av/QAV/D_ALLP
        for i in `seq 0 $Nmom_eq`; do imom=`printf %03d $i`; mkdir -p -v $dir_av/QAV/p$imom ; done

	args="-L=$L -T=$T  -dir=$dir_av -Nclust=${Nclust}"
	echo $args
	pwd
	#./D_pmAVERAGE $args -dir_plus=$dir_plus -dir_minus=$dir_minus
	#./Subst_Oa2g2_Z4F $args -NOa2=${NOa2}  -Uplaq=${Uplaq} -beta=${beta}  -gauge_fixing=${gauge_fixing} -gauge_action=${gauge_action} -Nf=${Nf}
	#for iOa2 in `seq 0 $NOa2m1`; do  for imom in `seq 0 $Nmom_eq`; do ./Calc_Z_oneP $args -imom=$imom -iOa2=$iOa2 ; done; done
	#./Z4FRGI $args -NOa2=${NOa2} -Nf=${Nf} -ainv=${ainv}
done #eList
fi

############################
## SEA CHIRAL LIMIT
###########################
if [ $DO_CHIRAL_SEA = "yes" ]
then
sh compile.sh write_allp_dir

#add the suffix of the rgi you want to analyse
#rgi=""
rgi="_RGI"

	# input sea
        echo $L > $dirinput/input_sea.dat
        echo $T >> $dirinput/input_sea.dat
        echo $Nclust >> $dirinput/input_sea.dat
        echo "$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/" >> $dirinput/input_sea.dat
	echo $rgi >>$dirinput/input_sea.dat

	for i in `seq 0 $Nmom_eq`; do imom=`printf %03d $i`; mkdir -p  -v $dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/a/QAV/p$imom/ ;  done
	mkdir -p  -v $dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/a/QAV/Z_ALLP/


	sea_massList=(0.0046 0.0078)
	sea_ensembleList=(3 5)
	# Mpcac
	#if [ $pm = "m" ] ; then sea_mpcacList=(-0.00585 -0.00823); fi
	#if [ $pm = "p" ] ; then sea_mpcacList=(0.00559 0.00829) ; fi
	# mPOL=sqrt(mu^2 +(za*mpcac)^2) mPOL=(mPOL,p + mPOL,m )/2	
	sea_mPOLList=(0.0064483 0.010167)
        seaNmass=${#sea_massList[*]}

        echo $seaNmass > $dirinput/mPOL_sea.dat
        for im in `seq 1 $seaNmass`; do ensemble=${sea_ensembleList[$(echo " $im -1 " |bc -l )]}; echo $ensemble >> $dirinput/mPOL_sea.dat; done
	for im in `seq 1 $seaNmPOL`; do  mPOL=${sea_mPOLList[$(echo " $im -1 " |bc -l )]}; echo $mPOL >> $dirinput/mPOL_sea.dat; done
	
#	time for imom in `seq 0 $Nmom_eq`;do for iOa2 in `seq 0 $NOa2m1`; do  sed "s/__IMOM__/$imom/g;s/__IOA2__/${iOa2}/g" ../PROGRAMS/Z4F_CHIRALSEA_template > Z4F_CHIRALSEA ; root -q -b -l Z4F_CHIRALSEA ; done	; done
	 dir_a="$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/a/"
	 args="-L=$L -T=$T  -dir=$dir_a  -rgiprint=$rgi -Nclust=${Nclust}  -NOa2=${NOa2}"
	./write_allp_dir $args  -
fi

#############################
## p2 fit & pm average
#############################
if [ $DO_P2FIT = "yes" ]
then

#rgi=""
rgiList=("_RGI")

# p2 cuts
p2min_lin=1.5
p2max_lin=2.2

p2min_av=1.1
p2max_av=1.23

cut=0.28

# create directory
for j in 1 2 ; do mkdir -p $dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/a/QAV/p2FIT/M$j ; done

for rgi in ${rgiList[*]}
do

        # input sea
        echo $L > $dirinput/input_fit.dat
        echo $T >> $dirinput/input_fit.dat
	echo "$p2min_lin $p2max_lin" >> $dirinput/input_fit.dat
	echo "$p2min_av $p2max_av" >>$dirinput/input_fit.dat
	echo "$cut" >>$dirinput/input_fit.dat
        echo $Nclust >> $dirinput/input_fit.dat
        echo "$dirdata/DATA/OUTPUT/${action}Nf${Nf}/b$beta100/Z4F/a/QAV" >> $dirinput/input_fit.dat
        echo $rgi >>$dirinput/input_fit.dat


	time for iOa2 in `seq 0 0`; do sed "s/__IOA2__/${iOa2}/g" ../PROGRAMS/Z4F_p2FIT_template > Z4F_p2FIT ; root -q -b -l Z4F_p2FIT ; done

done # rgi
fi

cd $dir
